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摘 要 : 球面 距离 ( 角 间 距 ) 计算 是 天 文 或 地 理学 中 极 常 用 的 计算 之 一 , 也 是 目标 查找 、 
锥 形 检 索 、 交 叉 证 认 等 方法 的 基础 。 数 学 上 ， 通 过 球面 几何 可 以 直接 计算 出 两 点 的 距离 ， 前 
人 也 已 经 推导 出 了 多 个 复杂 程度 不 一 进行 计算 。 但 是 由 于 计算 机 的 精度 有 限 , 在 进行 数值 计 
算 时 有 全 入 误差 ,导致 公式 计算 结果 出 现 偏差 ,对 几 个 常用 的 球面 距离 计算 公式 进行 了 考察 ， 
测试 并 对 比 它们 在 不 同 计算 环境 下 的 精度 与 优 缺点 。 此 外 还 展示 并 比较 了 几 种 常用 天 文 软件 
包 、 数 据 库 的 球面 距离 计算 方法 , 以 期 有 助 于 天 文 工 作者 选择 适合 自己 当前 需要 的 计算 方法 。 
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地 计算 出 来 。 由 于 
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的 长 度 , 通过 二 维 直 角 坐 标 使 用 勾 股 定理 可 以 很 容易 


只 涉及 乘法 、 加 法 与 开平 方 运算 ， 在 计算 机 中 也 可 以 保留 很 高 的 精度 。 而 


球面 距离 计算 则 要 复杂 得 多 。 球面 在 直角 坐标 系 中 实际 上 是 三 维 结构 ,除非 距离 相对 于 半径 


非常 之 小 ， 否 则 不 能 当 作 平面 处 理 。 


两 点 球面 距离 计算 是 天 文 或 地 理学 中 最 常用 的 计算 之 一 ， 是 目标 查找 、 锥 形 检索 、 交 又 


证 认 等 计 和 
算 一 遍 距离 ， 取 出 
球面 上 查找 与 茶点 的 
有 类 似 功能 "。 


中 最 基础 的 运 入 


。 比 如 查找 距离 某 个 位 置 最 近 的 山峰 , 需要 先 对 附近 的 山峰 都 计 


索取 得 与 其 最 近 
是 球面 距离 计算 。 
两 点 球面 距离 计算 ， 计 入 


点 之 间 较 小 的 那 段 圆 弧 的 长 度 ， 天 文 上 


Cangular separation) 或 角 吕 


E 离 最 小 的 一 个 。 而 锥 形 检 索 是 虚拟 天 文 台 的 数据 检索 协议 之 一 ,用 于 在 
E 离 小 于 指定 半径 的 目标 列表 , 几乎 所 有 提供 了 星 表 查 询 功 能 的 系统 都 
基于 位 置 的 星 表 交 叉 证 认 则 可 视 为 是 对 一 个 星 表 A 中 所 有 的 源 做 一 遍 锥 形 检 
的 另 一 星 表 的 源 ， 这 样 就 完成 了 A、B 两 个 星 表 的 交叉 ~“。 这 些 应 用 的 基础 都 


的 是 两 点 在 球面 上 的 最 短 距 离 , 即 球 心 与 两 点 组 成 的 大 圆 上 两 


习惯 使 用 该 短 圆 弧 对 应 的 大 圆 圆心 角 ， 称 为 角 间 距 


E 离 。 球 面 距离 需要 使 用 球面 几何 进行 求解 ， 要 大 量 使 用 三 角 


函数 。 而 计算 机 的 数值 计算 精度 有 限 , 在 对 三 角 函 数 、 反 三 角 函 数 进行 计算 时 容易 出 现 舍 入 
误差 。 多 个 函数 的 误差 积累 之 后 ,将 导致 严重 的 结果 偏离 ， 甚 至 无 法 得 到 正确 的 结果 。 前 人 


却 鲜 少 有 人 讨论 。 因 


已 经 推导 了 很 多 数学 公式 求解 球面 距离 , 但 在 何 种 情况 下 应 使 用 哪个 公式 , 公式 的 精度 如 何 ， 
此 ， 本 文 考察 天 文 技术 中 常用 的 球面 距离 计算 公式 ， 给 出 它们 的 算法 ， 


并 对 比 它们 的 精度 , 为 球面 距离 算法 的 选择 提供 参考 。 
1 球面 几何 方法 
球面 距离 计算 可 直接 通过 球面 几何 推导 出 来 。 较 常用 的 计算 公式 有 大 圆 (Great-circle) 


公式 、Haversine 公 式 等 。 另 外 还 有 Vincenty 公 式 被 应 用 到 Astropy 等 代码 库 中 ,但 其 复杂 
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度 过 高 ， 本 文 在 第 3. 1 节 中 简要 提 及 。 


设 有 两 点 赤道 坐标 为 (01; A) P2(02.02) ， 求 它们 的 球面 角 距离 d。 其 中 大 圆 公式 为 
JH OD ; 当 两 点 距离 很 小 时 , 也 有 人 使 用 简化 的 (2) 式 et, Haversine 公 式 如 (3) R. 


d = arccos [sin 9, - sin 02 + cosó, - cos ôz - cos (o4 — o3)] 


RENE 
d= le: — 03)- cos (21222) 


] cain . 2,01 — 024 a an 。 2 
d = 2 . arcsin 4/ sin“ (———— ) + cos ô - cosós - sin^( 


t (à) — 93)? 


Ot — €» 
2 


) (3) 


计算 机 中 通常 使 用 IEEE754 二 进 制 浮 点 数 算术 标准 处 理 浮 点 数 , 其 中 单 精度 32 位 (float 
或 single) 真 值 可 精确 到 小 数 点 后 6 9 位 , 双 精 度 64 位 (double ) 真 值 可 精确 到 小 数 点 后 15 17 
位 `。 本 文 使 用 C++ 语言 编写 程序 对 3 个 公式 的 精度 进行 了 对 比 ， 结 果 如 表 1。 本 文 使 用 的 单位 
ABE CO ， 表 中 有 两 行 数据 时 ， 上 方 的 值 为 单 精度 计算 的 结果 ， 下 方 值 为 双 精 度 计算 的 


R 


H, 


ZNR o 
R1 使 用 不 同 的 点 测试 常用 球面 距离 公式 的 精度 〈 单 位 : 度 ) 
Table 1 Accuracy testing at different points for widely used formulas (unit: 
degree) 

pl p2 简化 的 大 圆 公 式 大 圆 公 式 Haversine 

0. 00999999884516001 0. 00000000000000000 0. 00999999884516001 
0,0 0.01, 0 

0. 01000000000000000 0. 00999999998265327 0. 01000000000000000 

0. 00999999884516001 0. 00000000000000000 0. 00999999884516001 
0,0 0, 0.01 

0. 01000000000000000 0. 00999999998265327 0. 01000000000000000 

0. 00998573563992977 0. 00000000000000000 0. 00998573563992977 
180, 0 180. 01, 0 

0. 00999999999999938 0. 00999999998265327 0. 00999999999999938 

0. 00731059815734625 0. 00000000000000000 0. 00731059815734625 
42, 43 42. 01, 43 

0. 00731353701619125 0. 00731353703719986 0. 00731353701187370 

0. 00999939627945423 0. 00000000000000000 0. 00999939627945423 
42, 43 42, 43. 01 

0. 00999999999999302 0. 01000000001909975 0. 00999999999999302 

0. 01863496564328671 0. 00000000000000000 0. 01000121701508760 
0, 90 180, 89. 99 

0. 01862095887437574 0. 00999999998265327 0. 01000000000000639 
0,0 0.1,0 0. 09999999403953552 0. 10087054967880249 0. 09999999403953552 


? https://en.wikipedia.org/wiki/Haversine formula 
3 https://en.wikipedia.org/wiki/Vincenty9627s formulae 
^ https://en.wikipedia.org/wiki/Single-precision floating-point format 
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0. 10000000000000001 0. 10000000000019954 0. 10000000000000001 


由 表 1 可 见 ， 单 精度 下 , WALERAN, KAANE EE, BNT CERNI 
入 误差 , 而 距离 较 大 时 误差 也 较 大 。 双 精度 下 ， 大 圆 公 式 能 够 返回 结果 , 但 比 Haversine 公 式 
的 精度 要 差 。 大 圆 公式 的 简化 形式 多 数 时 候 与 Haversine 公 式 的 精度 相当 。 这 3 个 公式 共有 的 
一 个 问题 ， 它 们 均 有 可 能 略微 超过 准确 值 ， 当 严格 限定 一 个 边界 值 的 时 候 ， 可 能 导致 漏 源 。 
因而 在 实际 应 用 时 ， 需 要 考虑 将 边界 值 上 略微 外 扩 ， 以 将 因为 计算 精度 而 漏 掉 的 源 包含 在 内 。 

值得 注意 的 是 , 大 圆 公式 的 简化 形式 在 极点 附近 时 误差 非常 大 。 虽 然 它 在 上 述 公 式 里 是 
计算 复杂 度 最 低 的 ， 并 且 在 两 极 之 外 、 距 离 较 小 时 精度 与 Haversine 相 近 ， 但 是 在 实际 使 用 
该 公式 时 仍 应 当 非 常 小 心 。 


2 直角 坐标 系 方法 

在 基于 亦 道 坐标 系 进行 求解 之 外 , 也 可 以 使 用 直角 坐标 系 , 通过 三 维 坐 标 向 量 对 两 点 球 
面 距离 进行 计算 。 可 视 天 球 半径 为 单位 1， 以 天 球 球 心 为 三 维 直 角 坐 标 系 原点 ， 赤 道 坐 标 系 
北 天 极点 方向 为 z 轴 方向 , 天 赤道 面 上 赤 经 为 0 点 的 方向 为 x 正 向 , 遵循 右手 坐标 系 与 二 2 相交 
的 轴 为 z 匡 方 向 。 从 而 可 将 赤 经 、 杰 纬 转换 为 三 维 直 角 坐 标 (% y, 力 ， 其 转换 关系 如 (4-6) 


XU, 
£ = cos(6).cos(a) (4) 
y = cos(ô)-sin(a) (5) 
= sin(o) (6) 


这 样 两 点 距离 计算 可 以 直接 使 用 (8) 式 进行 计算 。 该 公式 可 简单 地 从 Haversine 公 式 中 
推 得 ， 其 中 Haversine 公 式 根 号 中 的 部 分 (7 式 ) 实际 上 是 两 点 的 三 维 空间 直线 距离 的 一 半 的 


ô; — ô: 9,01 — Q: 
A = si? ( : 7 2 ) + cos 61 - cos ó; - sin? ( 一 一) 

1 _ 

= ga [(£2 — 21)? + (ya — u1)? + (22 — 21)°] (7) 
2 二 
= sin (5) 

JL 1 - - - _ 
C d = 2.arcsinv A = 2- arcsin zi V (za — t1)? + (yo — 1)? + (22 — 4)? (8) 


直角 坐标 系 的 计算 精度 方面 ， 本 文 依旧 使 用 C++ 语 言 在 单 精度 及 双 精 度 两 个 环境 下 对 公 
式 进 行 计算 ， 并 与 Haversine 的 结果 对 比 ， 结 果 如 表 2。 从 表 2 可 以 看 出 ， 直 角 坐 标 系 方法 的 
精度 与 Haversine 几 乎 完全 一 致 ， 而 直角 坐标 系 方 法 有 时 甚至 优 于 Haversine 〈 见 42, 43 两 行 
数据 ) 。 


表 2 Haversine 与 直角 坐标 系 计 算 结 果 对 比 (Hw: E) 


Tab. 2 Result comparison between Haversine and Cartesian method (unit: degree) 


pl p2 Haversine 直角 坐标 系 方法 

0. 00999999884516001 0. 00999999884516001 
0, 0 0. 01, 0 

0. 01000000000000000 0. 01000000000000000 


0, 0 0, 0. 01 0. 00999999884516001 0. 00999999884516001 
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0. 01000000000000000 0. 01000000000000000 

0. 00998573563992977 0. 00998573563992977 
180,0 180. 01,0 

0. 00999999999999938 0. 00999999999999938 

0. 00731059815734625 0. 00731242587789893 
42, 43 42. 01, 43 

0. 00731353701187370 0. 00731353701187507 

0. 00999939627945423 0. 00999987311661243 
42, 43 42, 43. 01 

0. 00999999999999302 0. 00999999999999879 

0. 01000121701508760 0. 01000121701508760 
0, 90 180, 89. 99 

0. 01000000000000639 0. 01000000000000639 

0. 09999999403953552 0. 09999999403953552 
0,0 0.1,0 

0. 10000000000000001 0. 10000000000000001 


直角 坐标 系 方法 在 进行 距离 比较 时 还 可 以 进一步 进行 优化 ， 如 CO) 式 : 设 两 点 直角 坐 
PEAP Gn isi), pa(rs. ya 22)。 其 优点 是 可 以 预先 计算 x、y、z 与 threshold ( 设 最 大 角 
ERDO, 从 而 避免 在 实际 进行 距离 计算 时 使 用 非常 耗 时 的 三 角 函数 , 而 只 需要 三 次 减法 、 
三 次 乘法 、 两 次 加 法 、 一 次 比较 , 大 大 减少 了 计算 量 , 非常 适用 于 大 规模 数据 计算 。 相应 的 ， 
它 的 缺点 是 需要 额外 占用 存储 空间 来 保存 x、y、z 三 个 值 ， 在 对 内 存 占用 要 求 较 苛刻 的 环境 
中 可 能 会 成 为 问题 。 因 而 ， 实 际 的 公式 选择 ， 仍 应 视 数据 规模 、 计 算 环 境 而 定 。 另 外 ， 基 于 
三 维 向 量 的 计算 方法 ， 还 有 法 线 向 量 的 形式 ， 但 其 计算 过 程 比 本 节 要 复杂 ， 本 文 将 在 第 3. 3 
节 中 进行 简要 描述 。 


SR 
22.4 = (z2-— 21)? + (ya — i) + (22 — a)? < 2? - sin? (7) = threshold 


3 计算 库 

许多 天 文 软件 包 也 包含 了 球面 距离 计算 功能 。 本 部 分 选择 了 几 个 安装 、 使 用 较 便捷 
用 也 比较 广泛 的 库 , 演示 了 其 使 用 方法 ， 并 进行 测试 。 由 于 难以 分 辨 软件 包 中 实际 使 用 的 
单 精度 还 是 双 精 度 , 此 部 分 将 不 再 对 此 进行 区 分 , 而 直接 与 Haversine 双 精度 结果 进行 比较 。 


mm E 


3.1 Astropy 

当前 天 文 界 已 经 非常 流行 使 用 Python 来 进行 数据 处 理 ， 其 中 Astropy” 为 Python 的 推广 
传播 功 不 可 没 。Astropy 进行 球面 距离 计算 需要 使 用 astropy.units 及 
astropy. coordinates. SkyCoord 计 算 。 调 用 方法 如 下 ， 其 中 纹 两 点 之 间 的 角 距 离 〈 单 位 是 
度 ) : 


& = SkyCoord (raikunits. rad, deci*units. rad) 


c = SkyCoord (raskunits. rad, decskunits. rad) 

d = cu. separation (c). deg 

为 公平 起 见 ， 在 64 位 而 ndows 10 中 的 64 位 python3. 6. 3 重新 实现 了 Haversine 算 法 ， 并 与 
Astropy 的 距离 计算 函数 进行 了 对 比 。 对 比 结果 见 表 3， 可 以 看 到 两 种 方法 基本 一 致 ， 不 分 作 
仲 。 


站 Í 


(9) 


表 3 Windows 64 位 Python3 环境 下 Haversine 与 Astropy 计算 结果 对 比 〈 单 位 : E) 


Tab. 3 Result comparison between Haversine and Astropy on Windows 64bit Python3 (unit: 


degree) 
pl P2 Haversine Astropy 
0,0 0.01,0 0.01 0.01 
0,0 0,0. 01 0.01 0.01 
180,0 80. 01,0 0. 009999999999999377 0. 009999999999990905 
42, 43 42. 01, 43 0. 007313537011873697 0. 007313537011872698 
42, 43 42, 43. 01 0. 009999999999993016 0. 009999999999994572 
0, 90 80, 89. 99 0. 010000000000006394 0. 010000000000006394 
0,0 0.1,0 0.1 0.1 


阅读 Astropy 的 源 代 码 可 发 现 其 距离 计算 函数 separation 实 际 调 用 的 是 
astropy/coordinates/angle utilities. py 中 的 angular separation(lonl, latl, lon2, 
lat2) 函数 ， 其 代码 实现 使 用 Vincenty 公 式 ， 即 C100 式 。 该 式 在 所 有 位 置 的 距离 计算 都 更 
稳定 , 但 计算 复杂 度 也 更 高 ”。 从 上 述 计算 对 比 结果 来 看 , Astropy 并 没有 处 处 都 比 Haversine 
的 精度 更 高 ， 如 (180, 0) 一 行 ; 当然 ，Astropy 在 (42, 43) 上 的 精度 也 比 Haversine 更 好 。 但 
是 两 者 的 差距 都 非常 小 ， 几 乎 可 以 忽略 不 计 。 从 比较 结果 上 看 ，Haversine 己 经 相当 精确 。 
当然 ， 知 是 对 计算 稳定 性 有 更 高 要 求 ， 可 以 考虑 Vincenty 公 式 。 


Y : 2 N 。 c aie e 2 
[cos 0? - sin (a4 一 ao) + [cos ô - sinó? — sin, - cos ôz - cos (o4 — 3)] 
d — arctan ————— - z (10) 
sin 04 - sin 3 + cos 04 - cos ô> - cos (a4 — Q2) 


3. 2 HTM 


DEZH Hierarchical Triangular Mesh, HTM) ""JA& —ROSEERTEL EU E JJ. 336 
归 三 角 分 割 方法 。 分 层 三 角 网 格 目前 公开 提供 下 载 的 软件 包 中 包含 了 C# 代 码 库 及 一 个 SQL 
Server 扩 Æ 包 > 其 中 Spherical.Htm.Sgl.cs 文件 中 包含 了 函数 
fDistanceEq(ral, decl, ra2, dec2) ， 可 用 于 计算 距离 。 需 要 注意 的 是 ，fDistanceEq 的 返回 
值 单位 是 arcmin， 不 是 degree。 针 对 分 层 三 角 网 格 ， 本 文 使 用 C# 实 现 了 Haversine 函 数 并 与 
SphericalHTM v3. 1. 2 的 fDistanceEq 函 数 进 行 对 比 。 如 表 4， 其 结果 基本 相同 。 查 看 源 代码 
发 现 fDistanceEq 实 际 上 也 使 用 了 Haversine 公 式 进 行 计算 , 只 是 多 了 degree 转 arcmin 的 转换 。 
除了 C# 程 序 中 可 以 直接 调用 fDistanceEq 外 ，HTM 软 件 包 提供 的 SQL Server 数 据 库 扩展 
包 中 也 带 有 此 函数 ， 在 数据 库 中 使 用 非常 方便 。 


5 基于 Astropy 2.0.2， 参 考 网 址 
https://github.com/astropy/astropy/blob/master/astropy/coordinates/angle utilities.py 
6 http://www.skyserver.org/htm/index.html 


表 4 CHEMBET Haversine E HTM fDistanceEq 计算 结果 对 比 (单位 : E) 


Table 4 Result comparison between Haversine and HTM fDistanceEq on C# (unit: degree) 


pl p2 Haversine HTM fDistanceEq 
0,0 0.01,0 0.01 0.01 

0,0 0, 0. 01 0.01 0.01 

180,0 80. 01,0 0. 00999999999999938 0. 00999999999999938 
42, 43 42. 01, 43 0. 0073135370118737 0. 00731353701187507 
42, 43 42, 43. 01 0. 00999999999999302 0. 00999999999999879 
0, 90 80, 89. 99 0. 0100000000000064 0. 0100000000000064 
0,0 0.1,0 0.1 0.1 

3.3 SLALI1B 与 SOFA 


SLALIB 即 Starlink library of positional astronomy routines ， 使 用 Fortran 77 进 


行 实现 , 其 目标 是 让 天 文 工 作 者 更 简便 快捷 地 编写 精确 且 可 靠 的 天 体 测 量程 序 。 虽然 它 是 由 
Fortran 77 实 现 的 ， 但 是 已 经 有 人 为 其 编写 了 Python 接口 pyslalib ， 
Python 环境 中 获得 Fortran 77 的 计算 精度 。 


Pyslalib 中 的 slalib. sla dsep 直 接 对 应 了 SLALIB 中 
距离 计算 函数 ， 在 64 位 的 Ubuntu 16. 04. 4 操作 系统 中 使 用 64 位 Python3. 5. 2 环境 将 其 与 


Haversine 函 数 的 计算 结果 进行 了 对 比 ， 如 表 5。 


程序 计算 结果 与 Windows 一 致 ， 表 


因而 可 以 在 已 有 的 


的 sla DSEP (Al, B1, A2, B2) 球面 


首先 可 以 看 到 Linux 版 的 Python Harversine 


II 


明了 Python3 在 不 同 的 操作 系统 中 的 表现 一 致 。 其 次 


pyslalib 的 计算 结果 ,除了 (42, 43) 两 行 有 轻微 差别 , 与 Haversine 完 全 一 致 , 精确 也 非常 高 。 


查看 SLALIB 的 源 代码 ， 可 以 看 至 


与 第 2 节 不 同 , SLALIB 进 一 步 将 其 转换 为 法 线 向 量 , B 55722 , 再 使 | 


使 


士 FH 
结果 


天 体 测 量 中 常用 


o 


向 量 的 主要 优势 是 向 量 代 数 可 取代 部 分 三 角 函 数 的 计算 ， 
持 较 好 的 精度 ， 并 且 在 边界 点 〈 如 南北 天 极 、0-360 度 边界 ) 也 无 异常 
的 SOFA 库 ”的 球 
Standards of Fundamental Astronomy， 是 IAU 建 立 并 维护 的 一 个 权威 的 


计算 函 


1 其 使 用 了 第 2 节 


中 的 直角 坐标 系 三 维 向 


量 作为 计算 单元 。 但 


包含 有 ANSI C 及 Fotran 77 两 个 版 本 。 表 5 的 最 右 侧 列 出 了 双 精 度 位 C+ 


j (11) 式 ” 进行 计 算 。 
其 数值 计算 稳定 性 更 好 , 能 保 


o 


数 iauSeps 也 同样 使 


T4 x r5 
d = arctan ex Tal 
N1 T19 


7 http://stsdas.stsci.edu/cgi-bin/gethelp.cgi?slalib.sys 
8 https://github.com/scottransom/pyslalib 

9 https://en.wikipedia.org/wiki/N-vector 

10 http://www.iausofa.org 
11 http://www.iausofa.org/2018 0130 C/SeparationPA.html 


JY QD3xX. SOFAH 
基本 天 文学 标准 库 ， 
程序 调用 SOFA 计 算 的 


表 5 Linux 64 位 Python3 环境 下 Haversine 5E pyslalib 及 双 精 度 C++ 版 本 SOFA 计算 结果 


对 比 ( 单 


位 : 度 ) 


Table 5 Result comparison among Haversine and pyslalib on Linux 64bit Python3 and 


double precision C++ version SOFA (unit: degree) 
pl P2 Haversine pyslalib SOFA 
0,0 0. 01,0 0. 01 0. 01 0. 01000000000000000 
0,0 0, 0. 01 0. 01 0. 01 0. 01000000000000000 
180, 0 80. 01,0 0. 009999999999999377 0. 009999999999999377 0. 00999999999999938 
42, 43 42. 01, 43 0. 007313537011873697 0. 007313537011875421 0. 00731353701187542 
42, 43 42, 43. 01 0. 009999999999993016 0. 0099999999999971 0. 00999999999999710 
0, 90 80, 89. 99 0. 010000000000006394 0. 010000000000006394 0. 01000000000000639 
0,0 0.1,0 0.1 0.1 0. 10000000000000001 


4 数据 库 


关系 型 数据 库 均 支持 SQL 语句 ， 因 而 可 以 在 不 同 的 系统 中 直接 使 | 
于 每 种 数据 库 使 用 的 数学 函数 、 计 算 精 度 略 


式 进行 


E 离 计算 。 由 


SQL 语句 实现 Haversine 
有 差别 ， 本 文 在 Microsoft 


公 


SQL Server 2008、PostgreSQL 9.4.5、MySQL 5. 7. 18 上 进行 了 测试 。 测 试 结果 结果 如 表 6， 
从 表 6 可 以 看 到 ，3 个 数据 库 的 差距 不 大 ，MySQL 默 认 保留 的 精度 更 多 一 些 。 
表 6 不 同 数据 库 中 使 用 Haversine 公式 计算 距离 的 精度 对 比 〈 单 位 : 度 ) 


Table 6 Accuracy testing for pure SQL Haversine on different database (unit: degree) 


pl p2 MS SQL Server PostgreSQL MySQL 

0,0 0.01,0 0. 00999999999999926 0. 01 0.01 

0,0 0,0.01 0. 00999999999999926 0. 01 0. 01 

180, 0 80. 01,0 0. 00999999999999811 0. 00999999999999938 0. 009999999999999377 
42, 43 42. 01, 43 0. 00731353701187193 0. 0073135370118737 0. 007313537011873697 
42, 43 42, 43. 01 0. 00999999999999697 0. 00999999999999302 0. 009999999999993016 
0, 90 80, 89. 99 0. 0100000000000109 0. 0100000000000064 0. 010000000000006394 
0,0 0.1,0 0. 0999999999999995 0.1 0.1 


除了 Haversine 公 式 , 数据 库 中 使 用 直角 


标 系 方法 实际 上 更 


p 


方便 。 数 据 库 对 I/0 进 


T 


了 


大 量 优化 ， 可 以 更 好 地 对 数据 i 


4.1 PostgreSQL 数 据 库 插件 


数据 库 中 


减少 使 月 


插件 支持 球 


El cg 


专门 为 天 文 检索 设计 的 ， 
果 进 行 了 对 比 。 需 要 时 


行 计算 ， 


距离 计算 。 


除了 可 以 直接 使 用 公式 之 外 ， 也 
日 者 编程 的 麻烦 。 如 PostgreSQL 数 据 库 有 PostGIS”、pgsphere 、 
其 中 PostGIS 是 专 为 地 理 信息 系统 设计 的 ， 虽 然 也 可 以 在 天 文 上 使 
用 ， 但 是 需要 将 米 、 干 米 等 单位 转换 回 弧度 ， 太 过 繁琐 低 效 。 而 pgphere、Q3C” 、H3C 则 是 
更 为 简洁 。 下 文 演示 了 3 种 插件 的 使 用 方法 ， 表 7 对 这 3 个 插件 的 结 
昌 度 进行 计算 , 而 pgsphere 需 要 先 转 成 弧度 进 


FE 意 的 是 83C 和 H3C 直 接 使 月 


再 将 结果 转 回 角度 。 


SELECT a ral,b decl,c ra2, d dec2, 

q3c dist(a,b,c, d) q3c, 

h3c dist(a,b,c, d) h3c, 

DEGREES (spoint (RADIANS (a), RADIANS (b) ) €-^spoint (RADIANS (c), RADIANS (d))) pg 

FROM (VALUES (0.0, 0. 0, 0. 01, 0. 0) (0. 0, 0. 0, 0. 0, 0. 01), (180. 0, 0. 0, 180. 01, 0. 0), 
(42. 0, 43. 0, 42. 01, 43. 0), (42. 0, 43. 0, 42. 0, 43. 01), (0. 0, 90. 0, 180. 0, 89. 99), (0. 0,0. 0,0 


.1,0.0) AS v(a,b,c, d) 


可 以 使 用 


行 调度 。 因 而 ， 直 角 坐 标 系 方法 虽然 需要 增加 x、y、z 3 
个 字段 , 但 并 不 会 给 数据 库 增加 太 大 的 负担 , 还 可 以 大 大 降低 计算 的 复杂 度 ， 更 快 地 取得 海 


量 数据 的 计算 结果 。 


种 数据 库 扩展 插件 进行 相关 计算 ， 


Q3C". H3C^ 5E E 


$ 7 PostgreSQL 数据 库 Q3C. H3C. pgsphere 插件 的 球面 距离 计算 精度 比较 (单位 : RE) 


Table 7 Accuracy testing for different PostgreSQL extensions (unit: degree) 


pl p2 Q3C H3C pgsphere 

0.0,0.0 0. 01,0.0 0. 01 0. 01 0. 00999999998265327 
0.0,0.0 0. 0, 0. 01 0. 01 0. 01 0. 00999999998265327 
180. 0, 0. 0 180. 01, 0. 0 0. 00999999999999091 0. 00999999999999938 0. 00999999998265327 
42. 0, 43. 0 42. 01, 43. 0 0. 0073135370118727 0. 0073135370118737 0. 00731353703719986 
42. 0, 43. 0 42. 0, 43. 01 0. 00999999999999801 0. 00999999999999302 0. 0100000000190997 
0. 0, 90. 0 180. 0, 89. 99 0. 0100000000000064 0. 0100000000000064 0. 00999999998265327 
0.0,0.0 0.1,0.0 0.1 0.1 0. 1000000000002 


从 三 者 的 结果 来 看 , 精度 大 抵 相 当 , 没有 太 大 差距 , 在 极点 的 计算 也 都 没有 太 大 的 偏差 。 


经 阅读 源 代 码 ， 可 以 看 到 h3c_dist 函 数 使 用 了 Haversine 公 式 ; q3c dist K Zt si H 
Haversine 的 一 种 变 体 ; pgsphere 的 情况 则 较为 复杂 ， 


12 https://postgis.net/ 
13 https://github.com/akorotkov/pgsphere 

14 https://github.com/segasai/q3c 

15 http://cds.u-strasbg.fr/resources/doku.php?id-h3c 


局 在 距离 大 时 使 用 


w- 


T KE 公式 ， 距 


较 小 时 使 用 直角 坐标 公式 。 从 计算 结果 来 看 ， 这 3 个 插件 都 不 失 为 成 熟 可 用 的 扩展 ， 有 具体 要 
使 用 哪个 ， 还 要 看 没有 其 他 的 需求 。 如 Q3C、H3C 均 提供 了 join 函数 〈 分 别 是 q3c_join、 
h3c join) ， 对 两 星 表 交 叉 证 认 进行 了 优化 。 而 pgsphere 提 供 了 丰富 的 球面 几何 计算 ， 如 球 
面 形状 的 面积 计算 ， 球 面 形 状 的 交 、 并 计算 等 。 


5 总 结 


本 文 探讨 了 大 圆 公 式 、 简 化 的 大 圆 公 式 、Haversine 公 式 等 球面 距离 计算 方法 ， 并 对 比 
了 它们 在 计算 机 中 进行 单 精度 及 双 精 度数 值 计 算 的 结果 。 在 此 基础 上 , 引申 出 了 基于 三 维 直 
角 坐 标 系 的 距离 计算 方法 。 结 果 表 明 大 圆 公 式 在 单 精度 下 舍 入 误差 非常 严重 ， 而 大 圆 公 式 的 
简化 公式 在 两 极 或 两 点 距离 较 大 时 误差 很 大 。Haversine 公 式 与 直角 坐标 系 方法 的 精度 相近 、 
都 非常 高 ， 而 直角 坐标 系 方法 的 计算 量 较 小 ， 适 合 大 数据 的 计算 。 此 外 ， 还 需要 注意 所 有 这 
些 公式 在 计算 时 ， 有 时 候 计 算 结果 会 略微 大 于 准确 值 。 因 而 在 实际 应 用 时 ， 可 能 需要 将 边界 
值 取 得 略 大 一 点 ， 以 将 因为 舍 入 误差 漏 掉 的 源 包含 在 内 。 

已 有 软件 包 如 Astropy、HTM、SLALIB、SOFA 以 及 部 分 数据 库 扩 展 也 将 球面 距离 计算 模块 
包含 在 内 。 它 们 所 使 用 的 计算 方法 不 尽 相 同 ， 结 果 精 度 也 非常 高 。 当 需要 进行 球面 距离 计算 
时 ， 可 以 直接 使 用 这 些 库 ， 而 无 须 自行 实现 过 于 复杂 的 公式 ， 以 免 重 复工 作 且 易 因 考虑 不 周 
而 导致 出 错 。 
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